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Abstract 

The competition between reptation and Rouse Dynamics is incorporated in 
the Rubinstein-Duke model for polymer motion by extending it with sideways 
motions, which cross barriers and create or annihilate hernias. Using the Density- 
Matrix Renormalization-Group Method as solver of the Master Equation, the 
renewal time and the diffusion coefficient are calculated as function of the length 
of the chain and the strength of the sideways motion. These new types of moves 
have a strong and delicate influence on the asymptotic behavior of long polymers. 
The effects are analyzed as function of the chain length in terms of effective 
exponents and crossover scaling functions. 

PACS 61.25.Hk, 05.10.-a, 83.10.Kn 

1 Introduction 

A dilute solution of linear polymers in a gel provides the ideal case for reptation. The 
gel is a rigid network of obstacles which forces the polymer to find its way slithering 
through the maze. Effectively the polymer moves inside a tube of pores which changes 
only by growing and shrinking at the ends of the tube. It is not a great step to replace 
the network of the gel by a regular lattice. The regularity of the lattice still keeps 
the motion of the polymer random, because the ends randomly leave or enter cells 
of the lattice. A big step is the reduction of the motion to a stochastic process of 
hopping units. It certainly can not be justified on the level of monomers, because 
neighboring monomers are strongly correlated in their motion. For this purpose the 
notion of reptons has been introduced: blobs of monomers of the size of the correlation 
length [1]. Seeing the polymer as a sequence of reptons, permits to consider the units 
of motion as uncorrelated, with the only proviso that they do not separate too far in 
order to preserve the integrity of the polymer. 

Rubinstein [2] designed an elementary model for reptation, as a chain of slack and 
taut links connecting the reptons. A slack link describes two successive reptons in 
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the same cell and a taut link two in nearest neighbor cells. By allowing only these 
configurations of reptons and only moves between them, a simple model for reptation 
results. Duke [3j enriched the model by biasing the hopping of reptons by an external 
field, thus modeling the experimental situation of gel electrophoresis. Of course the 
model misses important aspects of polymer dynamics, such as the hydrodynamic inter- 
actions and even more importantly the requirement of self- avoidance, which infiuences 
the universal properties [HE]. One can incorporate self-avoidance in the model, which 
makes however the analysis an order of magnitude more difficult. An excuse for leaving 
out this aspect is that mutual exclusion is for reptons less severe than for monomers, 
because the reptons are loosely packed blobs of monomers which can interpenetrate 
each other. 

The eternal dilemma is to choose between being realistic and keeping it simple. 
The chemical culture opts for being realistic and deals with specific properties, the 
physical culture leans towards simplicity and aims at generic properties. It was de 
Gennes' contribution to polymer physics to show that properties of long polymers do 
not depend on the specific composition, in particular the dependence on the length of 
the polymers is governed by universal exponents. In this sense the Rubinstein-Duke 
(RD) model has shown to catch the essential physics of reptation in spite of the crude 
approximations made. 

We will approach the problem from the physical perspective and deal with the 
universal properties of the RD model, but investigate a richer class of motions than 
treated sofar. In the standard RD model only interchange of slack and taut links is 
permitted. It means that the length stored in a slack link moves in the direction of a 
taut link, interchanging the slack and taut link. A move that would perfectly fit in the 
spirit of the RD model, is the change of two consecutive slack links into two taut links. 
It corresponds to three reptons in the same cell of which the middle one escapes to a 
neighboring cell. Such a move does not cross a barrier. There was a practical reason 
to exclude this possibility, because it destroys the dimensional reduction, which was 
pointed out by Duke. From a physical point of view, the formation of "hernias" does 
not seem to infiuence the universal properties and in this paper we investigate this 
issue. Another optional move is the interchange of two taut links. As we will see this 
means the crossing of a barrier by the chain. If the barriers, posed by the obstacles were 
infinitely high, these processes would be strictly forbidden. But barriers are not perfect 
and therefore it is worthwile to investigate the infiuence of finite barriers. Moreover 
we will see that barrier crossing in the RD model does not lead to much change, but 
the combination with hernia creation and annihilation has the drastic effect of crossing 
over from reptation to Rouse dynamics. 

Since the standard RD model already does not permit an exact solution, we have 
to rely on numerical methods to analyze the extended model. The most common 
method is simulation of the system but this is less suited for our goal, because the 
crossover between the two types of dynamics occurs for rather long chains which are 
hard to simulate accurately. We will employ the technique of finite size analysis, which 
requires very accurate data to be successful. In this paper we use a method, based 
on the analogy between the Master Equation and the Schrodinger equation, by which 
the temporal evolution of the probability distribution of the chain configurations cor- 
responds to the evolution of the wave function. The Master operator corresponds to 
the Hamiltonian of a one- dimensional spin chain, for which the very efficient Density- 
Matrix Renormalization- Group Method (DMRG) has been designed by White pj. The 
model remains a one-dimensional quantum problem, irrespective the lattice in which it 
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is embedded, because the chain itself is a hnear structure. Apphcation of the DMRG 
method to the chain dynamics on the lattice is by now standard, but to perform success- 
fully calculations in a 3d embedding lattice requires an optimal use of the symmetries 
of the model in order to keep the basis set of states to a practical size. In the appendix 
we outline how we exploit the symmetries of the 3d lattice. 

In the next section we describe the extension of the RD model and the corresponding 
Master Equation. We focus on two properties: the renewal time r and the diffusion 
coefficient D and determine them directly from the Master Operator. The renewal 
time is the time needed for the chain to assume a new configuration, which has no 
memory of the original one. It is found from the gap in the spectrum of the Master 
Operator. The Master Equation always has a trivial eigenvalue 0, corresponding to the 
stationary state. Any other initial state ultimately decays towards the stationary state 
and the slowest relaxation time (the inverse of the gap) is the renewal time. The gap 
decays with a negative power z of the length N of the chain, such that r ~ N^. The 
zero field diffusion coefficient D is related to the drift velocity in a weak driving field 
and decays as a power A^~^. The approach to this asymptotic behavior is the main 
issue of this paper. 

Due to the new types of hopping, the dimensionality d of the embedding lattice plays 
a non-trivial role. We report, for the first time, on calculations in d = 3. They became 
possible through subtile use of the symmetries of the model, which are discussed in the 
appendix. The subsequent sections contain the results for the renewal exponent z and 
the crossover functions, which describe the data for all lengths N and strengths of the 
transition rates for barrier crossing and hernia creation/annihilation. The results for 
the diffusion coefficient D and its exponent x are calculated from a linearization of the 
Master Equation with respect to the driving field. The exponents z and x are linked 
through the mean square displacement of the wandering chains. In the discussion we 
comment on the results and explain why the crossover in gels is different from that in 
polymer melts. 

Paessens and Schiitz |9] have also extended the RD model by including "con- 
straint release" in the hopping rates. In our language this is a mix of hernia cre- 
ation/annihilation and barrier crossings. We comment their calculation in the discus- 
sion. 

Earlier [lOj we have performed a similar investigation for the cage model (in d = 2), 
with similar conclusions as the present study. Investigation of the RD model elucidates 
in how far the crossover is model independent. 

2 The model 

The model consists of a chain of + 1 reptons located in the cells of a (hyper)cubic 
lattice. They are connected by A^ links, labeled by Y = {yi, ■ ■ ■ , i/n)- The links take on 
the value i/i = (slack) or any of the 2d vectors which connect a cell to its neighbors 
(taut). The corners of the squares in c? = 2 or the edges of the cubes in d = 3 are 
barriers for the chain. The reptons can move in three ways: 

1. The standard RD moves. A repton between a slack and taut link moves to the 
neighboring cell thereby interchanging the slack and taut link. For these moves 
no barriers have to be overcome. A move is illustrated in Fig. [H The strength 
of the hopping rate for RD moves sets the time scale and is therefore put equal 
to 1. 
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Figure 1: Rubinstein-Duke moves. Their rate is set to 1. 



2. The barrier crossings, of which an example is shown in Fig. [2l This is an 
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Figure 2: Barrier crossings with hopping rate c 



interchange of two taut hnks connected to the same repton that jumps over the 
barrier. The strength of the transition rate for such a move is taken to be c. 

3. The creation of a hernia is the change from two consecutive slack links of which 
the middle repton jumps to a neighboring cell. The annihilation is the reverse 
process. An example is shown in Fig. [31 Hernia creation and annihilation go 
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Figure 3: Hernia creation a) to d) and annihilation d) to a). Both processes go with 
rate are h. 



with the rate h. 

The statistics of the model is governed by the Master Equation for the probability 
distribution P(Y,t), where Y stands for the complete configuration {yi, ■ ■ ■ .Vn)- It 
has the general form 



(9P(Y,t) 
dt 



^ [iy(Y|Y')P(Y', t) - iy(Y'|Y)P(Y, t)] = M(Y, Y')P{Y', t). (1) 



The W^s are the transition rates of the possible motions that we have indicated in the 
above list. The matrix M combines the gain terms (in the off-diagonal elements) and 
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the loss terms (on the diagonal). M is the sum of matrices, for each repton one 



M(Y,Y') = E^^(Y,YO, (2) 

i=0 

where the sum runs over the reptons starting with the tail repton i = to the head 
repton i = N. The internal reptons induce transitions between two configurations 
which differ in two consecutive yi, the external (head and tail) repton change only 
Un viz. yi. If we view the links as "bodies" the problem is equivalent with a one- 
dimensional many body system with two body interactions between nearest neighbors. 

The matrix M is asymmetric because the transition rates are biased by a factor 
B = exp(e/2) for the moves in the direction of the field and by for the reverse 
process, e is a dimensionless parameter representing the strength of the field. M is 
a stochastic matrix since the sum over each column vanishes. Therefore M has an 
eigenvalue and the corresponding right eigenvector is the probability density of the 
stationary state. All other eigenvalues are negative. The smallest in magnitude is the 
gap, giving the slowest decay to the stationary state and the inverse of the gap we take 
as the definition of the renewal time. The diffusion coefficient D is calculated from 
an infinitesimal small driving field. The field induces a drift Vd and via the Einstein 
relation 

the diffusion coefficient results. It is determined by expansion of the Master Equation 
in powers of e. 

M = Mo + eMi + ---, P(Y) = Fo(Y) + ePi(Y) + --- (4) 

which leads to the equations 

MoPo = 0, MoPi = -MiPo. (5) 

The first equation is trivially fulfilled by a constant Pq(Y), since the matrix Mq is 
symmetric and the right eigenvector becomes equal to the trivial left eigenvector. The 
right hand side of the second equation is a known function of the configuration. Thus it 
yields the linear perturbation Pi. The drift velocity is an average over the distribution 
[7]. The linear term in e of the drift velocity involves the terms Pq and Pi. With these 
terms we can calculate the linear term in Vd and find with ([3]) the diffusion coefficient 
D. 



3 The Exponent z for the Renewal Time 

The easiest way to obtain the exponent of the relation r = A^^ is to make a log-log plot 
and determine the slope. In a previous publication [8J it was shown that this is rather 
misleading for the present problem. A much more sensitive check is to compute local 
exponents zn according to 

_lnr(A^ + l)-lnr(A^-l)_ rflnr 

~ ln{N + 1) -\n{N -1) ~ ZbiV' ^ ' 
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Figure 4: The influence of the basis size m on the renewal time exponent for two 
combinations of parameters: c = 0, h = 0.1 and c = 0.1, h = 0. 



which gives 2; as a function of the chain length A^. The function is the basic 
ingredient for further analysis. 



First we check how accurate it can be obtained from a DMRG calculation of which 
the convergence is determined by the basis size m. In Fig. IHwe show examples of poor 
and excellent convergence. Poor convergence occurs for the combination c = 0, h = 0.1, 
where we have no barrier crossings but substantial hernia creation/annihilation. For 
a length = 20 the size of the basis, even as large as m=160, still has an influence 
and for longer chains this becomes stronger. Thus it is difficult to deduce from these 
data the effective exponent for chains longer than N = 20. In the second combination 
c = 0.1, h = the convergence is perfect for much smaller bases and for much longer 
chains, as the upper curve in the picture demonstrates. In fact the case with no 
barrier crossing at all is the only combination where convergence is a problem, as 
subsequent pictures will show. We blame the lack of convergence to the fact that 
hernia creation/annihilation without some barrier crossing leads to a large weight for 
configurations with many hernias and therefore to a short end-to-end distance of the 
chain. These are a-typical configurations and the DMRG procedure has difficulty to 
find an adequate basis to represent the gap state. 

Some other noticeable points are: 

1. The curves have still not reached the asymptotic value for the values of of the 
order of 100. Thus a log-log plot would suggest a higher value than the reptation 
value 3. This slow approach to the asymptotic value has been identified as the 
main reason for the discrepancy between the theoretical reptation exponent z=3 
and the measured higher values [8j. We get a better grip on the asymptotics 
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when we discuss the crossover. 



2. The curves do not indicate a tendency towards the Rouse exponent z = 2. This 
illustrates the point made in the introduction that the two mechanisms have to 
ai 




Figure 5: The renewal time exponent zj<i for /i = and a set of values c. 

The next Fig. [5] shows a set of curves for h = and a set of values c. Note that 
the effective exponent is quite sensitive to the value of c, but all curves do not show 
reptative behavior as was mentioned earlier for the dependence of a single point in 
parameter space. For larger values of c, a maximum in the effective exponent zn seems 
to develop for larger and larger A^. As the maximum can easily be interpreted as a 
saturated asymptotic value in a log- log plot, these corrections to scaling, still present 
for very long chains, are very important for assessing the correct asymptotic behavior. 
This feature makes it neccesary to do a finite size analysis in order to get a grip on the 
region in A^ where the behavior changes. 

The standard argument is that hernia creation/annihilation does not change the 
reptative character of the chain motion because they leave the backbone of the chain 
invariant, which is the collection of taut links after the chain has been successively 
stripped from its hernias. The backbone only changes by refreshment at the ends of 
the chain. Also barrier crossing seems to be, as a single mechanism, ineffective. It 
changes the backbone, but not the number of taut links in a certain direction, since 
in a barrier crossing taut links are only interchanged in position along the chain. We 
may call the properties of the chain, which do not change by internal motion, "quasi 
invariants". So one needs both a non-zero c and a non-zero h to remove these quasi 
invariants. The cooperation of hernia creation/annihilation and barrier crossing is a 
intricate mechanism. Therefore we concentrate first on the situation that one of them 
has a finite strength and the other becomes small. 
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The case where one of the two, c or h, is fixed at a finite value, gives a simple 
crossover from reptation to Rouse dynamics when the chain grows with A^. As example 
consider first a fixed value c = 0.1 and h varying and small, for which the local exponent 
zn is given in Fig. [6l The reverse situation is plotted in Fig. [7] with zn for h = 0.1 
and a set of values of c. The two figures are strikingly similar. For small values of the 
parameter c, (h) the chain seems to show reptative behavior but turns over towards 
the Rouse exponent z = 2 for longer chains. It is remarkable that even in Fig. [7] the 
values for very small c show this trend, while we know from the previous section that 
for c = the calculation is poorly convergent. 

Anticipating the asymptotic values of the two regimes, the following representation 
is adequate for the renewal time (for fixed h). 

T{N,c) = N'g{c'N). (7) 

The idea is that all curves of e.g. Fig. [7] are represented by a single curve g{x). Thus 
we have plotted in Fig. [8] the data for tN~^ as function of c^N for the fixed value 
h = 0.5 and varying c, with an assumed value 6 = 0.58. This exponent is determined 
by trial and error to get the maximum collapse of the data on a single curve. The figure 
shows indeed a nice data collapse but it hides a subtlety which we can uncover using 
the properties of the crossover function g{x). The function g{x) should be expandable 
for small arguments as 

g{x)=go + gix^ (8) 

and for large arguments as 
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Figure 7: The renewal time exponent zn for /i = 0.1 and a set of values c. 



Then r goes as A^^ for vanishing c and as A^^ for — oo (at non-zero c). Inserting 
the asymptotic behavior ([9]) into ([7]) we obtain 



where the dots refer to corrections of order In Fig. [9] we have made a plot of the 
limit of ln(r/A^^) vs Inc. The values of the vertical axis are extrapolated to iV — » oo, 
which corresponds to the first two terms of (JTOl) . We first check whether the basis of 
states is large enough, such that we have no systematic errors due to a too small basis. 
Then we inspect whether t/N"^ has a well defined limiting value for N oo. (The 
values should approach the limiting point in a fairly linear way.) If the curves in Fig. 
[9] had a straight slope, a well defined value of the crossover exponent 9 follows. As one 
observes, there rather is a constant slope for small values of c and another one for the 
larger values of c. Now crossover may only be expected in the limit of c — * 0, which 
gives a slope in the neighborhood of 6' ~ 0.5. We could therefore discard the behavior 
for larger c, as not being described by crossover, but this contrasts the findings for 
the cage model, where the crossover formula applies for practically the whole range of 
c. We show the data also for larger values of c because we find it intriguing that this 
region is also representable by a crossover function, albeit with a different crossover 
exponent. 



With one of the parameters h or c fixed it is the other parameter which controls 
the crossover. The real challenge is to find a representation where both mechanisms 
feature. We have not been able to find a simple expression which accurately accounts 
for arbitrary combinations of h and c. We gain some insight in the combined action of 
h and c by approaching the limit h = c = along a radial line h/c = r. For fixed r 



ln(r/A^2) = ln^_i-einc+---. 



(10) 



5 Crossover along lines h/c = r 
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Figure 8: The crossover scaling function g{x) for h = 0.5 and 6 = 0.58. 



we have again a single parameter which, in combination with A^, provides a crossover 
scaling variable, such that we can use the scenario of the previous section to analyze 
the data. In Fig. [10] we give the crossover exponent 6* as a function of In c, (In h) and 
for some values of the parameter ratio h/c. For very small values of Inc, (Ink) the 
slope of the line is compatible with the "universal" exponent 6 = 0.5. However for 
larger values another exponent seems to emerge of the order 6 = 0.85. Note that the 
curves in Fig. [10] run quite parallel, which means that r only enters in the offset given 
by g_i in (fTOj) . To show this point in more detail we have plotted in Fig. [11] the lines 
for larger values of h (or c). The window where the large exponent 6 = 0.85 applies is 
more than an order of magnitude for a fixed ratio h/c. 

6 The Diffusion Coefficient D 

The diffusion coefficient has been determined by the linearization (jSj), which gives the 
linear response of the drift velocity with respect to the driving force. We do not repeat 
the analysis for the diffusion exponent x, since diffusion and renewal time are closely 
related. If the center of the chain has drifted over a distance of the order of the end- 
to-end distance \/N, the chain has renewed itself. The mean square displacement due 
to diffusion during a renewal time equals Dr. So one has the connection 

Dt^N for N oo. (11) 

This relation implies for the exponents z — x = 1. We have tested this relation and in 
Fig. [12] we show the values of 2;7v ~ 1 and x^v for the same set of parameters, one for the 
small value h = c = 0.0001, where the behavior is more reptative and one for the larger 
value h = c = 0.1, where the exponents tend to Rouse dynamics. See also Fig. [13] 
which shows that the crossover exponent 6 for renewal and for diffusion are practically 
the same in the domain where it could be calculated with reasonable accuracy. 
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Figure 9: The crossover exponent 9 as deduced from fllOl) . It is given by the slope of 
the curve. 

7 Two-dimensional results 

We have shown the results for embedding dimension d = ?>. As mentioned the gen- 
eral expectation is that the embedding dimension has little influence on the universal 
properties. We are now in a position to verify this statement, since we have made 
extensive calculations both in ci = 2 and = 3. Indeed we come to the conclusion that 
the results agree qualitatively. To show an example we plot in Fig. [13] the d = 3 and 
d = 2 curves for 6 for h = c. The trends are the same, but the value of 6 in the "large" 
parameter regime is definitely larger for d = 3 than for d = 2. 

8 Discussion 

As a follow up of the study of crossover in the cage model, the extended RD model 
gives by and large the same picture: with growing length the chain crosses over from 
reptation to Rouse dynamics. The story is here more complicated because the two extra 
types of hopping, barrier crossing and hernia creation/annihilation, have to assist each 
other in order to get Rouse dynamics for long chains. This makes a comprehensive 
representation of the data in one scaling expression complicated. We have investigated 
the crossover behavior along lines in the c, h plane. 

The underlying idea of crossover is that there are two competing time scales. One 
is the diffusive time scale N'^/Dc, which is the time needed for a perturbation to 
diffuse along the chain inwards. We note that here not the overall diffusion D but the 
curvilinear diffusion coefficient applies, which decays as N~'^. This timescale leads 
to a renewal time ~ A^'^. The other time scale is the time needed to renew the chain by 
the combined action of hernia creation/annihilation and barrier crossing. If one of the 
parameters c or /i is large enough, the other is the limiting factor. If c is the smaller 
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Figure 10: The crossover exponent 9 for ratios of h/c. 



one, the timescale due to c equals N/c and if h sets the rate, it is N/h. The fastest of 
the two time scales sets the overall rate and therefore crossover occurs when they are 
equal i.e. when c ~ A^~^, or ~ A^~^, whichever is the smaller parameter. This leads 
to a crossover exponent 9 = 1/2. We see this trend in the numerical data, but the fact 
that we have to go to really small values of c (or h), and therefore to correspondingly 
large A^, prevents to see this "universal" crossover exponent in a clear way. 

On the other hand we observe for larger values of the parameters also crossover 
behavior, with different crossover exponents 9. This change in value could be a demon- 
stration of corrections to scaling, just as the exponents zjq or are rather far from 
their asymptotic values, also when crossover plays no role (as e.g. Fig. [5] shows). 
Clearly the renewal of the chain by sideways motion is still slow enough, even when c 
and h are of order unity, such that the competition with the diffusive renewal deter- 
mines the character of the dynamics. 

As pointed out earlier for the cage model [10], the crossover differs from the common 
scenario for polymer melts, where the crossover is in the opposite direction: from Rouse 
dynamics to reptation [11]. In the melt reptation results for the longer chains because 
the restriction in motion of the polymer, due to the presence of others, becomes more 
severe the longer the polymer is. We have argued that such crossover in melts can 
be understood from sideways motion that have a rate depending on the length of the 
chain. If the renewal time is taken as indicative for the lifetime of a barrier, the sideways 
motion would have a rate ~ N~^^ in the melt. The combined scaling parameter cN'^ 
then would shrink as grows. Since we always find zn > 2 reptation prevails for long 
chains in the melt. One would have to do a self consistent calculation, as carried out 
by Paessens and Schiitz, to make this argument quantitative [9]. 

Paessens and Schiitz [9] have also extended the RD model with rates that depend 
on the length of the chain. Their aim is to see the influence of "constraint release" on 
finite chains. The constraint release that they allow is in our language a mix of hernia 
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Figure 11: The crossover exponent 9 for ratios of h/c. 



creation and annihilation and barrier crossing. But not all types of barrier crossing 
that we allow, are permitted in their model. So, it is a bit difficult to make a clear 
comparison between their findings and ours. The interesting point of their calculation 
is the requirement of selfconsistency: the rates determine the renewal time and the 
renewal time in turn infiuences the rates. To carry out this program accurately within 
the DMRG method is one of the challenges for further research. 

Acknowledgment The authors have benefitted very much from many stimulating 
discussions with Gerard Barkema. A.D. thanks Wroclaw Centre for Networking and 
Computing for access for their computing facilities (grant No. 82). 

A The Symmetries of the Master Operator 

It is easy to set up a DMRG without paying attention to the symmetries of the prob- 
lem. Then in = 3 each link can be in 2(i + 1 = 7 states, leading to 7^ configurations 
for iV + 1 reptons. The possible symmetries in the problem will give an equal proba- 
bility to many configurations. If the symmetries are not explicitely acknowledged, the 
symmetries get lost when the choice of basis states does not conserve the symmetry. 
This means that states, which are equivalent by symmetry, have to be chosen simul- 
taneously. Thus either one has to include a large number of states, which leads to 
impractical calculations, or one has to keep track of the symmetry in each step of the 
method, which requires a substantial extra amount of careful programming. However, 
since we want to extract the utmost out of the data, we have no choice and must 
optimize the symmetry. 

As we deal with the fieldless gap for the renewal time and with the fieldless equation 
(l5l), we can employ, in principle, the full symmetry group of the cube, which has 
48 elements. If we were to do an exact calculation, we could apply all the relevant 
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Figure 12: Comparison of the exponents zn — ^ and xn- 



symmetry operations to the wavefunction and so reduce the number of components. 
But that is not the way DMRG works. The configuration space is sphtted into a tail 
and head part and the wavefunction is improved by an optimal choice of basis states 
in one part using the density matrix induced by the other part. In order to keep the 
symmetry in the wavefunction, the chosen states have to have the symmetry, which 
implies that the density matrix must have the symmetry. That in turn implies that, 
in each stage of the calculation, the wavefunction of the whole chain must have the 
desired symmetry. Thus we have to know how to combine symmetry of the parts in 
order to get the symmetry of the whole. This is similar to combining angular momenta 
of particles in atomic physics to get the angular momentum of the total wavefunction. 
The "good quantum numbers" derive from a set of commuting symmetry operators. 

We can find at most 3 commuting operators within the cubic group, with some 
freedom of choice. The simplest would be to look to the reflection symmetry of the 
coordinate axes. For each of the 3 operations the wavefunction can be even or odd, 
giving 8 sectors labeled by the parities. The parities qualify as good quantum numbers. 
The groundstate (stationary state) is located in the sector which is even under all three 
refiections. Each of the other sectors contains an excited state and the smallest (in 
magnitude) is the gap. The parities of the parts can be easily be combined to parities 
for the total since they simply multiply. We have implemented this scheme, but it does 
not lead to very accurate results, which we blame to the rather unbalanced occupation 
of the sectors, when the most probable states of the density matrix are chosen. 

The most successful use of the symmetry comes from another choice of commuting 
symmetry operations. We put the field in the direction of the body diagonal and 
consider rotations around this diagonal. We may rotate over the angles = 0, 27r/3 or 
47r/3, leaving the problem invariant. Under a rotation the wavefunction is multiplied 
with a phase factor exp(i0), which qualifies also as a good quantum number, leading 
to 3 sectors. Rotations commute with simultaneous inversion of the coordinate axes. 
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Figure 13: Comparison of the crossover exponent 9 ior d — 2 and d — 3, both for the 
gap and the diffusion coefficient. 



doubhng each of the 3 rotation sectors. As first 3 sectors we take those which are 
invariant under inversion and the next 3 are odd under inversion. The groundstate is 
in the first and the gap is the fourth sector (invariant under rotation and odd under 
inversion). The advantage of these good quantum numbers is that we can combine 
the quantum numbers of parts (by simple multiphcation) to the same set of quantum 
numbers for the combination. For instance a part in sector 2 (0 = 27r/3) and one in 
sector 3 (0 = 47r/3) lead to a combination with = 0, which is therefore in sector 1. 
This use of symmetry gives good results, but it is not yet optimal. 

A refinement could be made by considering the interchange of the x and the y axis. 
This turns the first and fourth sector into itself and transforms sector 2 and 3 as well as 
4 and 5 into each other. Although there is no good quantum number associated with 
this operation, we could use this symmetry by splitting sector 1 and 4 into an even 
and odd part under the interchange. This leads to 8 "channels" , which partly coincide 
with the previous sectors. Combining a part in sector (channel) 2 with a part in 3 does 
give an overall state in sector 1, but one has to take even and odd combinations to get 
them in the even and odd channel corresponding to sector 1. This gives a substantial 
amount of extra programming in order to keep properly track of the channels. But the 
effort is rewarded, as it improves the accuracy which is needed for the delicate cases 
of the parameter space. This choice of sectors (channels) is efficient because the states 
of the density matrix, which are chosen as having the largest eigenvalues, are more or 
less evenly distributed over the channels. 
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